* Figures B.10, B.11, B.3 panel c)
* 29 Oct 2021

clear all
set more off

grstyle init
grstyle set plain, nogrid horizontal
grstyle color background white

*Choose figure width:
local fig_w= 5
local res = `fig_w' * 600
	

* Main calibration
import excel "$output_s/w3_test_of_model.xlsx", sheet("Sheet1") firstrow clear 

su psi_h, d
replace psi_h = . if psi_h < r(p5) | psi_h > r(p95)
su psi_l, d
replace psi_l = . if psi_l < r(p5) | psi_l > r(p95)
su A, d
replace A = . if A < r(p5) | A > r(p95)

***************************************************************
** Figure B.10

local T = 1000*(exp(2.333) - 1)
set level 95
	*Colored version
	twoway (lpolyci A k_bl_post_trans if k_bl_post_trans <= `T', bw(250) lcolor(navy) acolor(navy*0.2)) (lpolyci A k_bl_post_trans if k_bl_post_trans > `T' & k_bl_post_trans < 13000, bw(250) lcolor(navy) acolor(navy*0.2)), xline(`T', lcolor(maroon)) xsca(r(7500 13000)) xla(8000(1000)13000) ytitle("A (excluding 5% outliers)") yla(#6) xtitle("Capital at baseline (after transfer) in BDT", margin(t+2)) legend(off) name(Lpoly_A, replace)
	graph save Lpoly_A "$output/Lpoly_A.gph", replace
	graph export "$output/Lpoly_A.png", as(png) replace

	twoway (lpolyci psi_l k_bl_post_trans if k_bl_post_trans <= `T', bw(250) lcolor(navy) acolor(navy*0.2)) (lpolyci psi_l k_bl_post_trans if k_bl_post_trans >= `T' & k_bl_post_trans < 13000, bw(250) lcolor(navy) acolor(navy*0.2)), xline(`T', lcolor(maroon)) xsca(r(7500 13000)) xla(8000(1000)13000) ytitle("psi_l (excluding 5% outliers)") yla(#6) xtitle("Capital at baseline (after transfer) in BDT", margin(t+2)) legend(off) name(Lpoly_psi_l, replace)
	graph save Lpoly_psi_l "$output/Lpoly_psi_l.gph", replace
	graph export "$output/Lpoly_psi_l.png", as(png) replace

	twoway (lpolyci psi_h k_bl_post_trans if k_bl_post_trans <= `T', bw(250) lcolor(navy) acolor(navy*0.2)) (lpolyci psi_h k_bl_post_trans if k_bl_post_trans >= `T' & k_bl_post_trans < 13000, bw(250) lcolor(navy) acolor(navy*0.2)), xline(`T', lcolor(maroon)) xsca(r(7500 13000)) xla(8000(1000)13000) ytitle("psi_h (excluding 5% outliers)") yla(#6) xtitle("Capital at baseline (after transfer) in BDT", margin(t+2)) legend(off) name(Lpoly_psi_h, replace)
	graph save Lpoly_psi_h "$output/Lpoly_psi_h.gph", replace
	graph export "$output/Lpoly_psi_h.png", as(png) replace


	*TIFF Greyscale version

	*Panel A
	twoway (lpolyci A k_bl_post_trans if k_bl_post_trans <= `T', bw(250) lcolor(gs5) acolor(gs10)) (lpolyci A k_bl_post_trans if k_bl_post_trans > `T' & k_bl_post_trans < 13000, bw(250) lcolor(gs5) acolor(gs10)), xline(`T', lcolor(black)) xsca(r(7500 13000)) xla(8000(1000)13000) ytitle("A (excluding 5% outliers)") yla(#6) xtitle("Capital at baseline (after transfer) in BDT", margin(t+2)) legend(off) name(Lpoly_A, replace) xsize(`fig_w')
	graph export "$output/TIFF/Lpoly_A.tif", as(tif) replace width(`res')
	
	*Panel B
	twoway (lpolyci psi_l k_bl_post_trans if k_bl_post_trans <= `T', bw(250) lcolor(gs5) acolor(gs10)) (lpolyci psi_l k_bl_post_trans if k_bl_post_trans > `T' & k_bl_post_trans < 13000, bw(250) lcolor(gs5) acolor(gs10)), xline(`T', lcolor(black)) xsca(r(7500 13000)) xla(8000(1000)13000) ytitle("psi_l (excluding 5% outliers)") yla(#6) xtitle("Capital at baseline (after transfer) in BDT", margin(t+2)) legend(off) name(Lpoly_psi_l, replace) xsize(`fig_w')
	graph export "$output/TIFF/Lpoly_psi_l.tif", as(tif) replace width(`res')
	
	*Panel C
	twoway (lpolyci psi_h k_bl_post_trans if k_bl_post_trans <= `T', bw(250) lcolor(gs5) acolor(gs10)) (lpolyci psi_h k_bl_post_trans if k_bl_post_trans > `T' & k_bl_post_trans < 13000, bw(250) lcolor(gs5) acolor(gs10)), xline(`T', lcolor(black)) xsca(r(7500 13000)) xla(8000(1000)13000) ytitle("psi_h (excluding 5% outliers)") yla(#6) xtitle("Capital at baseline (after transfer) in BDT", margin(t+2)) legend(off) name(Lpoly_psi_h, replace) xsize(`fig_w')
	graph export "$output/TIFF/Lpoly_psi_h.tif", as(tif) replace width(`res')

***********************************************************************************************
** Figure B.11

	*Colored version:
	twoway (hist psi_l, freq w(0.005) col(navy%40)) (hist psi_h, freq w(0.005) col(maroon%40)), xtitle("psi_h and psi_l (excluding 5% outliers)", margin(t+2)) legend(order(1 "psi_l" 2 "psi_h") pos(1) ring(0)) name(hist_psi, replace)
	graph save hist_psi "$output/hist_psi.gph", replace
	graph export "$output/hist_psi.png", as(png) replace

	*TIFF Greyscale version
	twoway (hist psi_l, freq w(0.005) col(gs5)) (hist psi_h, freq w(0.005) fcol(gs10%50) lcol(gs10) lpattern(dash)), xtitle("psi_h and psi_l (excluding 5% outliers)", margin(t+2)) legend(order(1 "psi_l" 2 "psi_h") pos(1) ring(0)) name(hist_psi, replace) xsize(`fig_w')
	graph export "$output/TIFF/hist_psi.tif", as(tif) replace width(`res')

***********************************************************************************************
** Figure B.3 panel c)

	*Colored version:
	twoway (hist A, freq w(50) col(maroon%40)), xtitle("A (excluding 5% outliers)", margin(t+2)) name(hist_A, replace)
	graph save hist_A "$output/hist_A.gph", replace
	graph export "$output/hist_A.png", as(png) replace

	*TIFF Greyscale version
	twoway (hist A, freq w(50) fcol(gs10) lcol(gs5)), xtitle("A (excluding 5% outliers)", margin(t+2)) name(hist_A, replace) xsize(`fig_w')
	graph export "$output/TIFF/hist_A.tif", as(tif) replace width(`res')


* Calibration using waves 2&3, 3&4, 4&5
forv j = 3/5 {
    import excel "$output_s/w`j'parameters.xlsx", sheet("Sheet1") firstrow clear 
	
	su psi_h, d
	replace psi_h = . if psi_h < r(p5) | psi_h > r(p95)
	su psi_l, d
	replace psi_l = . if psi_l < r(p5) | psi_l > r(p95)
	su A, d
	replace A = . if A < r(p5) | A > r(p95)
	
	twoway (hist psi_l, freq w(0.005) col(navy%40)) (hist psi_h, freq w(0.005) col(maroon%40)), xtitle("psi_h and psi_l (excluding 5% outliers)", margin(t+2)) legend(order(1 "psi_l" 2 "psi_h") pos(1) ring(0)) name(hist_psi, replace)
	graph save hist_psi "$output/hist_psi_w`j'.gph", replace
	graph export "$output/hist_psi_w`j'.png", as(png) replace

	twoway (hist A, freq w(50) col(maroon%40)), xtitle("A (excluding 5% outliers)", margin(t+2)) name(hist_A, replace)
	graph save hist_A "$output/hist_A_w`j'.gph", replace
	graph export "$output/hist_A_w`j'.png", as(png) replace
}


** All in one
import excel "$output_s/w3_test_of_model.xlsx", sheet("Sheet1") firstrow clear  
keep hhid psi_h psi_l A
gen wave = 2

su psi_h, d
replace psi_h = . if psi_h < r(p5) | psi_h > r(p95)
su psi_l, d
replace psi_l = . if psi_l < r(p5) | psi_l > r(p95)
su A, d
replace A = . if A < r(p5) | A > r(p95)

save "$output/parameters_all.dta", replace

forv j = 3/5 {
    import excel "$output_s/w`j'parameters.xlsx", sheet("Sheet1") firstrow clear 
	su psi_h, d
	replace psi_h = . if psi_h < r(p5) | psi_h > r(p95)
	su psi_l, d
	replace psi_l = . if psi_l < r(p5) | psi_l > r(p95)
	su A, d
	replace A = . if A < r(p5) | A > r(p95)
	gen wave = `j'
	append using "$output/parameters_all.dta"
	sleep 500
	save "$output/parameters_all.dta", replace
}


save "$output/parameters_all.dta", replace

tw (kdens A if wave == 2) (kdens A if wave == 3) (kdens A if wave == 4) (kdens A if wave == 5), ytitle("Density") xtitle("A (excluding 5% outliers)", margin(t+2)) name(kdens_A, replace) legend(order(1 "Waves 1 and 2" 2 "Waves 2 and 3" 3 "Waves 3 and 4" 4 "Waves 4 and 5"))
graph export "$output/kdens_A_all.png", replace
tw (kdens psi_h if wave == 2) (kdens psi_h if wave == 3) (kdens psi_h if wave == 4) (kdens psi_h if wave == 5), ytitle("Density") xtitle("psi_h (excluding 5% outliers)", margin(t+2)) name(kdens_psi_h, replace) legend(order(1 "Waves 1 and 2" 2 "Waves 2 and 3" 3 "Waves 3 and 4" 4 "Waves 4 and 5"))
graph export "$output/kdens_psi_h_all.png", replace
tw (kdens psi_l if wave == 2) (kdens psi_l if wave == 3) (kdens psi_l if wave == 4) (kdens psi_l if wave == 5), ytitle("Density") xtitle("psi_l (excluding 5% outliers)", margin(t+2)) name(kdens_psi_l, replace) legend(order(1 "Waves 1 and 2" 2 "Waves 2 and 3" 3 "Waves 3 and 4" 4 "Waves 4 and 5"))
graph export "$output/kdens_psi_l_all.png", replace
